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Abstract 


DRDC  Atlantic  recently  acquired  a  Reson  8125  bathymetric  sonar.  This  is  a  very 
high  resolution  multibeam  sonar  that  requires  accurate  measurement  of  the  instan¬ 
taneous  attitude  and  position  of  the  vessel  in  order  to  properly  determine  seabed 
bathymetry.  During  early  deployments  of  the  instrument,  roll  and  heave  contami¬ 
nation  of  the  bathymetric  data  resulted  from  less-than  optimal  configuration  of  the 
system  ancillary  sensors.  In  addition,  the  data  show  evidence  of  unresolved  motions, 
mainly  flex  in  the  sonar  mounting  pole. 

In  the  absence  of  commercial  software  tools  to  deal  with  this  problem,  a  method  has 
been  developed  for  removing  some  of  these  artifacts  from  the  data.  It  is  primarily 
intended  for  use  in  salvaging  useful  information  from  the  early  data  sets  now  that  the 
Reson  system  setup  has  subsequently  been  improved,  however  it  is  likely  that  pole 
flexure  will  continue  to  be  an  issue  when  deploying  this  sonar  on  CFAV  Quest.  The 
algorithm  has  been  implemented  in  Matlab  (The  Mathworks,  Inc.),  and  a  student, 
David  Cowan,  has  written  a  version  for  IDL  (Research  Systems,  Inc.). 


Resume 


RDDC  Atlantique  a  recemment  fait  l’acquisition  d’un  sonar  bathymetrique  Re¬ 
son  8125.  II  s’agit  d’un  sonar  multi-faisceaux  a  tres  haute  resolution  exigeant  la 
mesure  precise  de  l’assiette  et  de  la  position  instantanees  du  navire  pour  les  leves 
bathymetriques  des  fonds  marins.  Dans  les  premieres  applications  de  l’instrument, 
la  contamination  des  donnees  bathymetriques  par  roulis  et  levee  etait  due  a  la  config¬ 
uration  moins  qu’optimale  des  capteurs  auxiliaires  du  systeme.  De  plus,  les  donnees 
revelent  la  presence  de  mouvements  non  resolus,  principalement  le  flechissement  du 
poteau  servant  au  montage  du  sonar. 

En  l’absence  sur  le  marche  de  logiciels  permettant  de  resoudre  ce  probleme,  on 
a  mis  au  point  une  methode  en  vue  de  supprimer  certains  de  ces  artefacts  des 
donnees.  Cette  methode  vise  principalement  a  sauvegarder  l’information  utile  des 
premiers  ensembles  de  donnees,  maintenant  que  Installation  du  systeme  Reson  a 
ete  amelioree,  mais  il  est  probable  que  le  flechissement  du  poteau  continue  de  faire 
probleme  une  fois  que  ce  sonar  sera  monte  a  bord  du  NAFC  Quest.  L’algorithme  a 
ete  mis  en  application  chez  Matlab  (The  Mathworks,  Inc.),  et  un  etudiant,  David 
Cowan,  en  a  prepare  une  version  pour  IDL  (Research  Systems,  Inc.). 
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Executive  summary _ 

Introduction 

High  resolution  bathymetric  sonars  have  direct  application  in  the  area  of  mine  de¬ 
tection  by  high  frequency  sonars  and  for  other  naval  applications  which  require 
high  accuracy  and  resolution  in  imaging  of  the  seabed  and  underwater  objects. 
When  DRDC  Atlantic  first  acquired  a  Reson  8125  multibeam  bathymetric  sonar, 
there  were  some  configuration  issues  which  led  to  contamination  of  the  collected 
bathymetry  data  by  roll,  heave  and  other  ship  motion  artifacts.  A  method  has  been 
developed  that  removes  some  of  these  so  that  useful  information  can  be  derived  from 
those  early  data  sets.  DRDC  Atlantic  does  not  possess  commercial  software  that 
can  do  this,  though  such  tools  probably  are  available. 

Significance  of  Results 

After  applying  the  cleaning  method  to  data  sets  collected  on  the  Scotian  Shelf  (Q267) 
and  over  mine-shaped  targets  in  Esquimalt,  smaller  scale  features  are  resolvable  in 
the  bathymetry  data  that  would  otherwise  be  masked  by  the  artifacts. 

Future  Plans 

Most  of  the  deployment  issues  which  compromised  the  earlier  bathymetry  measure¬ 
ments  with  the  Reson  system  have  been  resolved,  however,  artifacts  due  to  unmea¬ 
sured  motions  (mainly  flex  in  the  mounting  pole)  may  continue  to  be  present.  This 
data  cleaning  method  is  applicable  to  the  components  of  pole  motion  that  appear 
as  roll  and  heave,  and  so  may  continue  to  be  useful  for  cleaning  data  collected  with 
the  multibeam  sonar  in  the  future. 


Crawford,  Anna  M.  2003.  Removing  Roll  and  Heave  Artifacts  from  High-Resolution 
Multibeam  Bathymetric  Data.  TM2003-243.  DRDC  Atlantic. 
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Introduction 

Les  sonars  bathymetriques  a  haute  resolution  ont  une  application  directe  dans  le 
domaine  de  la  detection  des  mines  par  des  sonars  a  hautes  frequences  ainsi  que  dans 
dautres  applications  navales  qui  requiarent  une  precision  et  une  resolution  extremes 
dans  limagerie  des  fonds  marins  et  des  objets  sous-marins.  Lorsque  RDDC  Atlan- 
tique  a  fait  l’acquisition  d’un  sonar  bathymetrique  multi-faisceaux  Reson  8125,  il  ex- 
istait  divers  problemes  de  configuration  occasionnant  la  contamination  des  donnees 
bathymetriques  par  des  artefacts  de  roulis,  de  levee  et  d’autres  mouvements  du 
navire.  On  a  mis  au  point  une  methode  en  vue  de  supprimer  certains  de  ces  artefacts 
pour  que  l’information  utile  puisse  etre  extraite  de  ces  ensembles  de  donnees.  RDDC 
Atlantique  ne  possede  pas  de  logiciel  disponible  sur  le  marche  qui  puisse  le  faire, 
bien  qu’il  en  existe  probablement. 

Signification  des  resultats 

Apres  L  application  de  la  methode  d’epuration  aux  ensembles  de  donnees  recueillis 
sur  le  Plateau  neo-ecossais  (Q267)  et  sur  des  cibles  en  forme  de  mines  a  Esquimalt,  il 
est  possible  d’effectuer  la  resolution  d’elements  a  echelle  inferieure  dans  les  donnees 
bathymetriques,  qui  seraient  autrement  masques  par  les  artefacts. 

Plans  futurs 

La  plupart  des  problemes  d’application  qui  compromettaient  les  premieres  mesures 
bathymetriques  effect uees  au  moyen  du  systeme  Reson  ont  ete  resolus,  mais  il  se 
pourrait  que  des  artefacts  dus  a  des  mouvements  non  mesures  (principalement  le 
flechissement  du  poteau  de  montage)  demeurent  toujours  presents.  Cette  technique 
depuration  des  donnees  s’applique  aux  composantes  de  mouvement  du  poteau  qui 
apparaissent  comme  du  roulis  et  de  la  levee  et  pourrait  done  s’averer  utile  pour 
l’epuration  des  donnees  recueillies  dans  l’avenir  au  moyen  du  sonar  multi-faisceaux. 


Crawford,  Anna  M.  2003.  Suppression  des  artefacts  de  roulis  et  de  levee  des  donnees 
bathymetriques  multi-faisceaux  a  haute  resolution.  TM2003-243.  RDDC  Atlan¬ 
tique. 
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1  Introduction 


A  fully  operational  multibeam  sonar  system  incorporates  a  large  number  of  com¬ 
plicated  sub-systems  and  ancillary  sensors,  all  of  which  need  to  be  understood  and 
properly  set  up  for  accuracy  of  the  resulting  bathymetry  measurements.  DRDC 
Atlantic  recently  acquired  a  Reson  8125  multibeam  sonar  which  is  a  fine  example  of 
such  a  complicated  system.  The  ancillary  sensor  suite  consists  of  a  heading  sensor 
(gyro),  a  Motion  Reference  Unit  (MRU),  a  sound  velocity  probe  (SVP)  and  a  Dif¬ 
ferential  Global  Positioning  System  receiver  (DGPS).  During  the  first  few  months  of 
operations  with  this  new  sonar,  some  of  these  sensors  were  not  properly  configured 
or  were  not  communicating  with  each  other.  The  result  of  this  was  contamination 
of  the  collected  bathymetry  data  with  artifacts  of  the  ship  motion  (mostly  roll  and 
heave).  In  addition,  the  data  shows  evidence  of  motion  of  the  sonar  head  relative 
to  the  ship  due  to  flex  in  the  sonar  mounting  pole,  which  is  not  measured  by  the 
motion  sensor. 

Over  the  course  of  several  operations  subsequent  to  the  initial  deployments  of  the 
system,  problems  with  the  configuration  have  been  eliminated  one  by  one  so  that 
currently  the  system  is  approaching  an  optimal  operating  state.  A  remaining  factor 
which  cannot  be  easily  remedied  is  the  flex  of  the  mounting  pole. 

Though  the  future  of  operations  with  this  high-resolution  instrument  promise  ex¬ 
citing  science,  we  are  left  with  survey  data  from  the  early  history  of  its  use  that 
contain  motion  contamination.  DRDC  Atlantic  does  not  currently  have  commer¬ 
cial  software  tools  that  can  deal  with  this  problem,  however  a  method  has  been 
devised  for  removing  the  roll  and  heave  components  of  these  artifacts.  This  was 
developed  and  written  in  Matlab  (The  Mathworks,  Inc.),  and  has  also  subsequently 
been  implemented  by  a  student,  David  Cowan,  in  IDL  (Research  Systems,  Inc.). 
This  method  has  allowed  extraction  of  some  interesting  smaller  scale  bathymetric 
features  from  the  motion-contaminated  data  sets. 
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2  Multibeam  Geometry  Considerations 


A  multibeam  sonar  system  measures  bathymetry  by  sensing  the  range  to  the  seabed 
along  each  beam,  then  knowing  the  angles  of  the  beams  and  the  attitude  and  po¬ 
sition  of  the  sonar  head,  determines  the  positions  of  the  beam  footprints  on  the 
seabed.  The  Reson  8125  has  240  half-degree  beams  arranged  in  a  fan  shaped  wedge, 
with  the  outermost  beams  at  ±60°  from  vertical  to  either  side  of  the  sonar  nadir 
(athwartships)  and  a  1°  fore-aft  beamwidth.  As  the  vessel  travels  forward,  a  swath 
along  the  seafloor  is  mapped  which  has  width  determined  by  the  water  depth.  For 
each  sonar  ping,  determining  where  each  of  the  240  beams  intersect  seabed  requires 
accurate  knowledge  of  the  attitude,  position  and  heading  of  the  sonar  head  at  the 
time  of  that  ping.  Furthermore,  accuracy  is  required  in  measurement  of  the  relative 
positions  of  the  sonar  head,  the  ship  motion  sensor,  and  the  center  of  gravity  (CoG) 
of  the  vessel  since  these  are  not  co-located.  Rotations  and  translations  of  the  vessel 
are  assumed  to  act  about  the  CoG,  but  are  measured  at  the  locations  of  the  MRU 
and  GPS  receiver  antenna,  and  the  resulting  motion  of  the  sonar  head  is  then  calcu¬ 
lated  using  geometry.  Figure  1  illustrates  a  3-dimensional  view  of  a  vessel  carrying 
a  multibeam  system,  including  the  sonar  head,  some  ancillary  sensors  and  the  CoG 
of  the  vessel. 

For  future  reference,  a  line  refers  to  the  set  of  seabed  positions  derived  from  one 
transmitted  ping  over  all  beams  (240  of  them  for  the  Reson  system).  Each  line  then 
contains  a  bathymetric  profile  generally  running  perpendicular  to  the  alongtrack 
direction  (the  ship  heading).  A  series  of  lines  is  compiled  as  a  swath. 

There  are  a  number  of  ways  that  inaccuracy  can  introduce  motion  artifacts  or  other 
errors  into  bathymetry  data. 

•  Inaccurate  CoG  or  Sensor  Offsets:  If  any  of  the  relative  positions  of  the 
sonar  head,  CoG,  GPS  receiver  and  MRU  are  incorrectly  determined,  then  the 
geometric  calculation  of  the  motion  and  position  of  the  sonar  head  from  the 
measured  attitude  and  position  will  give  incorrect  results. 

•  Sensor  Latency:  There  is  a  time  delay  between  measurement  of  ship  attitude 
by  the  MRU,  or  of  ship  position  by  the  GPS,  and  reporting  of  this  information 
to  the  sonar  processing  software.  This  can  result  in  previous  motion  or  position 
measurements  being  applied  to  the  current  sonar  ping  ranges  while  the  attitude 
and/or  position  has  likely  changed  in  the  meantime. 

•  Miscommunication  Between  Sensors:  In  order  to  operate  optimally,  some 
system  sensors  require  aiding  information  from  others.  The  MRU  in  particular 
needs  to  be  in  communication  with  the  heading  sensor  (gyro)  so  that  centrifugal 
acceleration  is  not  confused  with  other  motions. 

•  Sound  Speed  Profile:  The  sonar  head  has  a  co-located  sound  speed  probe, 
the  output  of  which  aids  in  beamforming.  If  the  sound  speed  in  the  surrounding 
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Figure  1:  3-dimensional  view  of  a  vessel-mounted  multibeam  system.  Rotation  and 
translation  of  the  vessel  are  about  the  CoG.  Position  is  measured  at  the  GPS  receiver 
and  attitude ,  at  the  MRU. 

water  column  differs  from  that  at  the  point  of  measurement  or  from  that  assumed 
by  the  processing  software,  the  difference  between  real  and  calculated  acoustic 
refraction  results  in  erroneous  positioning  of  the  beam  footprints  on  the  seabed. 

•  Unresolved  Motions:  If  there  is  relative  motion  between  the  attitude  or  po¬ 
sition  sensors  and  the  sonar  head,  this  is  not  measured  and  therefore  will  not  be 
included  in  the  beam  footprint  position  calculations.  As  previously  mentioned, 
an  example  of  this  is  flex  in  the  sonar  mounting  pole. 

There  are  other  more  subtle  causes  which  tie  into  one  or  more  of  the  above  factors. 
For  example,  the  motion  sensor  uses  Kalman  filtering  to  condition  the  attitude  out¬ 
put.  In  addition  to  needing  the  aiding  gyro  information,  it  also  requires  filter  settings 
that  are  appropriate  to  the  operating  environment.  Motion  in  Very  long  swell  waves 
may  not  be  properly  resolved  when  operating  with  the  settings  appropriate  for  short 
sea  waves. 

In  more  recent  deployments,  the  latency  and  sensor  miscommunication  issues  have 
been  resolved.  The  sonar  processor  can  utilize  a  Pulse  Per  Second  (PPS)  time  tag 
output  from  the  GPS  to  attach  the  correct  attitude  and  position  measurements  to 
the  appropriate  set  of  ping  beam  footprints.  This  eliminates  the  latency  problem. 
In  addition,  the  serial  communications  between  the  gyro,  motion  sensor  and  sonar 
processor  have  now  been  correctly  configured. 


DRDC  Atlantic  TM2003-243 


3 


Figure  2:  An  example  of  roll  and  heave  contaminated  bathymetry  data .  The  mean 
water  depth  is  68  meters. 

The  artifacts  evident  in  the  bathymetry  data  collected  during  the  first  few  deploy¬ 
ments  of  the  DRDC  Atlantic  Reson  8125  are  mainly  due  to  a  combination  of  three 
of  the  above-listed  factors  —  sensor  latency,  miscommunication  between  sensors 
and  unresolved  motion.  These  show  primarily  as  what  will  be  referred  to  as  roll 
and  heave  artifacts,  though  would  more  properly  be  called  roll-like  and  heave-like . 
Roll  artifacts  are  characterized  by  vertical  deviations  that  increase  in  amplitude  to¬ 
ward  the  outer  beams.  These  are  equivalent  to  rotation  of  the  ship  about  a  fore-aft 
axis,  with  deviation  of  opposite  sign  on  opposite  sides  of  center.  Heave  artifacts 
are  equivalent  to  vertical  motions  of  the  ship  and  have  a  constant  value  along  lines. 
Both  types  of  artifacts  can  be  recognized  by  the  alongtrack  oscillations  they  induce. 

Figure  2  illustrates  an  example  of  motion-contaminated  data,  mostly  suffering  from 
roll  artifacts.  The  outer  beams  show  deviations  of  up  to  a  meter  in  68  meter  water 
depth.  There  is  a  lesser  degree  of  heave  contamination  in  this  particular  example. 
This  data  was  collected  during  cruise  Q267  aboard  CFAV  Quest  on  the  Scotian  Shelf 
while  waves  were  running  3-4  meters  high.  Beam  footprint  positions  that  have  been 
rejected  by  the  sonar  processor  automatic  quality  control  algorithms  leave  holes 
(white  spaces)  in  the  surface  in  the  image. 
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3  The  Roll  and  Heave  Cleaning  Method _ 

The  algorithm  for  removing  roll  and  heave  contamination  proceeds  as  follows: 

•  The  sonar  processing  software  exports  the  calculated  positions  of  the  beam  foot¬ 
prints  (X,  Y,  Z  values),  along  with  line  and  beam  numbers,  in  a  formatted  ASCII 
file.  The  position  of  the  sonar  head  and  attitude  information  that  accompany 
each  ping  is  output  in  a  header  with  each  line  of  beam  footprint  coordinates. 

•  The  processing  script  reads  the  lines  of  beam  position  data  and  compiles  series 
of  slopes  and  intercepts  of  least-squares  linear  regressions  to  each  line,  relative 
to  the  sonar  head  position. 

•  The  series  of  line  slopes  and  intercepts  are  filtered  to  remove  the  high  frequency 
alongtrack  variation,  leaving  the  lower  frequency  (longer  length)  changes. 

•  For  each  line,  the  slope  is  adjusted  to  the  filtered  slope  by  rotating  about  the 
sonar  head  position,  then  the  intercept  is  adjusted  to  the  filtered  intercept  value 
by  shifting  the  whole  line  vertically.  All  three  coordinates  of  the  beam  footprint 
positions  for  a  line  are  altered  by  the  slope  adjustment. 

•  The  adjusted  beam  footprint  positions  are  output  as  a  new  ASCII  file  with  the 
same  format  as  the  original  file. 

A  Matlab  implementation  of  this  algorithm  is  presented  in  Annex  A.  It  is  likely  that 
there  are  commercially  available  software  packages  that  use  algorithms  equivalent 
to  this.  The  Ocean  Mapping  Group  at  University  of  New  Brunswick  (John  Hughes 
Clark),  for  example,  has  developed  tools  that  work  this  way. 

Figure  3  illustrates  the  series  of  line  slopes  and  intercepts  that  accompany  the  con¬ 
taminated  data  shown  in  Figure  2.  The  filtered  values  are  also  shown  (in  red).  The 
adjustments  to  the  intercept  are  minor,  at  most  about  15  cm.  The  longer  wavelength 
trends  in  slope  and  intercept  are  retained,  since  they  may  represent  real  bathymetric 
features.  The  small  scale  acrosstrack  variations  in  each  line  are  also  retained. 
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Figure  3:  Raw  and  filtered  line  slopes  and  intercepts  for  the  sample  of  data  shown  in 
Figure  2. 

4  Samples  of  Cleaned  Bathymetric  Data 

Figure  4  illustrates  the  output  of  the  cleaning  algorithm  applied  to  the  data  that 
was  shown  in  Figure  2.  Both  the  original  and  cleaned  data  are  also  shown  in  plan 
view  in  Figure  5. 

An  interesting  small  scale  feature  in  the  cleaned  data  shown  in  Figure  5  are  the 
bedforms  in  the  lower  elevation  (blue)  area  to  the  right  in  the  images.  These  are 
ripple  bedforms  in  the  coarse  sand  and  gravel  sediments  found  there,  with  wave¬ 
lengths  about  2  m  and  heights  about  0. 1-0.2  m.  The  boundary  between  the  higher 
(yellow-green)  and  lower  elevation  area  is  a  steep  scarp  with  slope  about  30°. 

Another  example  of  cleaned  bathymetry  data  is  shown  in  Figure  6.  The  upper 
panel  shows  cleaned  and  geo-referenced  bathymetry  data  (on  a  75  cm  by  75  cm 
grid),  while  the  lower  panel  shows  geo-referenced  sidescan  sonar  data  over  the  same 
area  (20  cm  by  20  cm  grid).  The  scour  holes  surrounding  the  mine-shaped  targets 
are  seen  in  the  bathymetry  data  in  a  mean  water  depth  of  about  60  meters.  A 
vertical  cylinder  target  just  above  the  center  of  the  images  is  itself  resolved  in  the 
bathymetry  data. 
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Figure  4:  The  results  of  cleaning  the  data  that  was  shown  in  Figure  2. 
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Alongtrack  line  number 


Uncorrected  bathymetry 


Roll  and  heave  artifacts  removed 


Acrosstrack  beam  number 


Figure  5:  Plan  views  of  the  contaminated  and  cleaned  bathymetry  data  shown  in 
Figures  2  and  4. 
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Figure  6:  A  cluster  of  mine-shaped  targets  in  cleaned,  geo-referenced  bathymetry 
data  and  in  sidescan  sonar  data  over  the  same  area. 
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5  User  Beware 


Removing  the  high  frequency  roll  and  heave  components  from  a  data  set  may  not 
resolve  all  the  problems.  Figures  7  and  8  illustrate  an  example  of  data  containing 
other  artifacts  that  cannot  be  removed  by  this  method.  The  cleaned  data  shown 
in  Figure  8  still  shows  oscillation  of  the  very  outer  beams  after  the  high  frequency 
variation  in  line  slope  has  been  removed.  The  residual  oscillation  has  the  same 
frequency  as  the  slope  variation  that  was  removed,  but  is  out  of  phase  by  90°.  Note 
that  the  oscillation  of  the  roll  artifact  which  has  been  removed  has  a  different  peak 
frequency  than  the  vessel  motions  (roll,  pitch  or  heave).  The  residual  may  be  due 
to  refraction  effects,  as  indicated  by  more  strongly  affecting  the  outer  beams. 

The  algorithm  as  presented  (see  Annex  A)  may  need  to  be  adjusted  to  suit  the 
particular  situation.  At  the  least,  the  filtering  of  the  slopes  and  intercepts  should 
be  tuned  to  the  motions  evident  in  the  data.  As  currently  implemented,  the  filters 
have  been  specified  with  window  lengths  defined  by  a  number  of  lines,  rather  than 
on  a  time  base.  This  assumes  that  the  sonar  ping  rate  and  the  frequencies  of  the 
unwanted  motion  contamination  do  not  change  significantly  over  the  duration  of  the 
set  of  lines  being  cleaned.  In  the  case  that  there  are  significant  measured  (real)  pitch 
or  roll  offset  values,  the  horizontal  position  of  the  sonar  head  may  fall  some  distance 
away  from  the  center  of  the  line  of  beam  positions  on  the  seabed.  Particularly  for  a 
large  pitch  offset,  correction  by  rotation  in  the  vertical  plane  containing  the  line  will 
result  in  distortion  of  the  adjusted  line,  especially  near  nadir.  Rotation  about  the 
sonar  head  position  is  in  fact  only  an  approximation  of  the  motion  being  removed 
from  the  data.  It  would  be  very  difficult  to  determine  the  true  center  of  rotation 
from  the  available  information. 

Finally,  this  method  is  applied  to  data  that  is  in  the  form  of  X,  Y,  Z  values,  referenced 
by  line  and  beam  numbers,  i.e.  data  that  is  in  a  very  raw  state.  This  technique 
cannot  be  applied  to  bathymetric  data  that  has  been  processed  such  as  by  gridding 
where  line  and  beam  numbers  (and  even  ship  track  and  heading)  are  lost. 
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elevation  elevation 


Figure  7:  Raw  Reson  bathymetry  data  showing  motion  contamination. 


Figure  8:  Results  of  cleaning  the  data  shown  in  Figure  7.  Note  the  persistence  of 
alongtrack  devation  in  the  outer  beams. 
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6  Conclusions 


The  method  that  has  been  presented  here  can  be  used  to  remove  undesirable  mo¬ 
tion  (roll  and  heave)  contamination  from  raw  multibeam  bathymetric  data.  Some 
examples  have  been  presented  showing  the  effectiveness  of  the  technique,  as  well  as 
a  few  warnings  to  potential  users.  The  algorithm  itself  is  presented  in  the  Annex. 

Application  of  this  technique  has  allowed  some  interesting  small  features  to  be  ex¬ 
tracted  from  bathymetry  data  that  would  otherwise  be  overshadowed  by  the  motion 
contamination.  The  Reson  8125  sonar  has  demonstrated  the  ability  to  resolve  rip¬ 
ple  bedforms  of  less  than  half  a  meter  height  in  almost  70  m  water  depth  in  rough 
seas.  The  example  of  the  mine-like  targets  shows  promise  for  the  system  in  future 
minehunting  research. 
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Annex  A:  Matlab  Implementation  of  the 
Algorithm _ 


The  Matlab  script  presented  below  is  not  the  most  computationally  efficient  method 
for  performing  these  calculations,  but  is  effective  and  illustrates  the  idea. 

Starting  at  the  point  where  the  data  has  already  been  imported  into  Matlab  vari¬ 
ables,  the  input  arrays  are  E,  N,  and  Z,  the  3-dimensional  (Easting,  Northing,  and 
depth)  beam  footprint  positions  with  indices  [line  number,  beam  number  from  1 
to  240],  and  vectors  are  refposE,  refposN  and  refposZ,  the  positions  of  the  sonar 
head  for  each  line.  In  the  arrays  E,  N  and  Z,  beam  positions  that  have  been  rejected 
by  the  sonar  processor  quality  control  algorithms  are  filled  with  NaN  (not  a  number). 

First,  the  series  of  line  slopes  and  intercepts,  relative  to  the  sonar  head,  are  de¬ 
termined.  R  are  the  ranges  of  the  beam  footprints  in  a  horizontal  plane  from  the 
horizontal  position  of  the  sonar  head.  Here,  R  is  made  negative  to  one  side  by 
flagging  negative  angles,  which  is  appropriate  for  headings  roughly  East-West  or 
West-East.  In  the  case  of  a  North-South/South-North  heading,  R  would  be  made 
negative  for  th  greater  than  +90°  or  less  than  -90°.  The  slope  and  intercept  of  a  line 
are  then  the  slope  and  intercept  of  the  set  of  values  (R,Z-refposZ) .  The  function 
find_slope  (presented  at  the  end  of  the  Annex)  calculates  a  linear  least-squares 
regression,  returning  the  slope,  intercept,  and  regression  coefficients.  Transforms  be¬ 
tween  cartesian  and  polar  co-ordinates  are  made  using  the  built-in  Matlab  functions 
cart2pol  and  po!2cart. 


for  il  =  l:size(E,l)  %  step  through  the  lines 

idvalid  =  f  ind(~isnan(Z  (il ,  :  )  )  )  ;  °/0  find  valid  values  (not  NaN) 

°/0  in  plan  view,  convert  to  polar,  with  origin  refposE, refposN 
[th,R] =cart2pol (E ( i 1 , idvalid) -refposE (il) ,N(il , idvalid) -refposN (il) ) 
flips  =  f ind(th<0) ;  °/0  make  R  negative  for  negative  angles 

R(flips)  =  -R(flips) ; 

°/0  now  find  slope  of  R,Z  (Z  is  always  negative) 

[slp(il)  ,intrcpt(il)  ,r(il)  ,s(il)]  =  ... 

find_slope(R,Z(il, idvalid)  -  refposZ (il) ) ; 

end; 


The  series  of  line  slopes  and  intercepts  are  then  filtered.  In  this  case,  a  5th-order 
Butterworth  filter  has  been  used  (built-in  function  butter)  that  has  a  width  of  20 
samples.  The  width  will  require  adjustment  according  to  the  motions  found  in  the 
data.  The  filter  is  applied  using  the  built-in  filtfilt  command,  which  results 
in  a  forward-backward  filter  operation  eliminating  phase  distortion  of  the  results. 
Before  filtering,  a  20-sample-long  pad  is  temporarily  added  to  the  start  and  end  of 
each  series  to  minimize  end  effects. 
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°/0  5th  order  butterworth  filter,  with  window  20  samples  wide 
[B,A]  =  butter (5, 0.05) ;  °/0  this  may  require  tweeking 

°/0  filter  intercepts,  filtfilt  is  a  zero-lag  filter 
fintrcpt  =  f iltf ilt (B, A,  .  .  .  °/0  pad  start  and  end  with  mean  values 
[ones (1 , 20) *mean(intrcpt)  intrcpt  ones (1 , 20) *mean(intrcpt)] ) ; 
fintrcpt  =  fintrcpt  (21 :  end-20) ;  °/0  now  remove  padding 
°/0  same  for  slopes 

fslp  =  f iltf ilt (B, A, [ones(l ,20)*mean(slp)  sip  ones (1 , 20) *mean(slp)] ) ; 
fslp  =  fslp(21 : end-20) ; 


Finally,  the  slope  and  intercept  of  each  line  are  adjusted  to  the  filtered  values. 


for  il  =  l:size(E,l)  °/0  step  through  lines 

idvalid  =  f ind(~isnan(Z(il ,  : ) ) ) ;  °/0  find  valid  values 

°/0  in  plan  view,  convert  to  polar,  with  origin  refposE,refposN 
[thl ,R] =cart2pol (E ( i 1 , idvalid) -refposE(il) ,N(il , idvalid) -refposN (il) ) ; 
flips  =  f ind(thl<0) ;  °/0  make  R  negative  for  negative  angles 

R(flips)  =  -R(flips) ; 

°/0  convert  to  polar  in  a  vertical  plane,  Z  rel.  to  sonar  head 
[th2,R2]  =  cart 2pol (are ,Z(il , idvalid) -ref posZ (il) ) ; 

°/0  then  convert  back  after  rotating  to  the  filtered  slope 
[R3,Znew(il, idvalid)]  =  pol2cart (th2-(atan (sip (il) -fslp (il) ) ) ,R2) ; 

°/0  and  calculate  new  E  and  N  along  original  acrosstrack  direction 
[Enew(il , idvalid) ,Nnew(il, idvalid)]  =  pol2cart (thl , abs (R3) ) ; 

°/0  now  make  E  and  N  relative  to  world  coordinates  again 
Enew(il , idvalid)  =  Enew(il , idvalid)  +  refposE(il); 

Nnew(il , idvalid)  =  Nnew(il , idvalid)  +  refposN(il); 

°/0  adjust  Z  co-ordinate  to  filtered  intercept,  restore  refposZ 
Znew(il , idvalid)  =  Znew(il , idvalid)  +  refposZ(il)  ... 

-  (intrcpt (il)-fintrcpt (il)) ; 

end; 

The  function  f  ind_slope: 


function  [slope , intercept ,r, s]  =  f ind.slope (X, Z) ; 

°/0  calculate  slope  and  intercept  from  linear  least  squares  regression. 

1  see  CRC  Math  Handbook  pgs  497-498  (27th  Ed.,  1985) 

°/0  input  -  X,  Z  =  same-sized  input  vectors  (your  data) 

°/0  output  -  slope,  intercept  =  of  linear  least  sq.  regression  to  your  data 
°/0  r,s  =  regression  statistics  (population  correlation  coeff 

°/0  and  standard  error) 
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°/0  trap  improperly  sized  input 

if  (length (X)  ~=  length (Z)),  error (’ Improper  input  vectors’);  end; 
N  =  length (X); 

slope  =  (N*sum(X.*Z)  “  sum(X) *sum(Z) ) / (N*sum(X . ~2)  -  sum(X)~2); 
intercept  =  mean(Z)  -  slope*mean(X) ; 
r  =  (N*sum(X.*Z)  -  sum(X) *sum(Z) ) / . . . 

sqrt ((N*sum(X. "2)  -  sum(X) "2) * (N*sum(Z . "2)  -  sum(Z)~2)); 
s  =  sqrt(sum((Z  -  (intercept+slope*X) ) . "2) / (N-2) ) ; 

return; 
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